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Abstract 

Recent experiments on coating flows and liquid drop impact both demonstrate that wetting failures 
caused by air entrainment can be suppressed by reducing the ambient gas pressure. Here, it is shown that 
non-equilibrium effects in the gas can account for this behaviour, with ambient pressure reductions increas¬ 
ing the gas’ mean free path and hence the Knudsen number Kn. These effects first manifest themselves 
through Maxwell slip at the gas’ boundaries so that for sufficiently small Kn they can be incorporated into 
a continuum model for dynamic wetting flows. The resulting mathematical model contains flow structures 
on the nano-, micro- and milli-metre scales and is implemented into a computational platform developed 
specifically for such multiscale phenomena. The coating flow geometry is used to show that for a fixed 
gas-liquid-solid system (a) the increased Maxwell slip at reduced pressures can substantially delay air en¬ 
trainment, i.e. increase the ‘maximum speed of wetting’, (b) unbounded maximum speeds are obtained as 
the pressure is reduced only when slip at the gas-liquid interface is allowed for and (c) the observed be¬ 
haviour can be rationalised by studying the dynamics of the gas film in front of the moving contact line. A 
direct comparison to experimental results obtained in the dip-coating process shows that the model recovers 
most trends but does not accurately predict some of the high viscosity data at reduced pressures. This 
discrepancy occurs because the gas flow enters the ‘transition regime’, so that more complex descriptions 
of its non-equilibrium nature are required. Finally, by collapsing onto a master curve experimental data 
obtained for drop impact in a reduced pressure gas, it is shown that the same physical mechanisms are also 
likely to govern splash suppression phenomena. 


1 Introduction 


The behaviour of liquid-solid-gas contact lines, at which an ambient gas is displaced by a liquid spreading 
over a solid, is the key element of many technological processes, ranging from coating flows used to apply thin 


liquid films to substrates (Weinstein & Ruschak 2004) through to 3D printers being developed to economically 


produce bespoke geometrically-complex structures (Derby 20101. Often, limitations on the operating devices are 


imposed by ‘wetting failure’ which occurs when the liquid can no longer completely displace the gas at the contact 
line. In coating technologies, this leads to the entrainment of gas bubbles into the liquid film which destroy 
its carefully tuned material properties whilst in drop-based flows entrainment of air under the moving contact 
line promotes the formation of splashes which, again, are usually undesirable. Consequently, understanding the 
physical mechanisms which control wetting failure and, where possible, delaying this occurrence by manipulating 
the material parameters, flow geometry, etc, to ever increasing wetting speeds is one of the key aims of research 
into dynamic wetting. 

Recent experimental studies have discovered new ways in which to delay the onset of wetting failures. 


Specifically, in both coating flows (Benkreira & Khan 2008 Benkreira & Ikin 2010) and drop impact phenomena 


(Xu et al. 2005 Xu, 20071 it has been shown that a sufficient reduction in the pressure of the ambient gas. 


around a factor of ten in the former case and even less in the latter, can suppress wetting failures. This result 


was unexpected, as until Xu et al. { 2005) it was assumed that the gas’ dynamics can either be entirely neglected 


or, as found in coating research, only influences the system’s behaviour through viscous forces which act in 


the thin gas film in front of the moving contact line (Kistler 1993). However, the significance of the viscous 
forces in the gas is usually characterised by the viscosity ratio between the gas and the liquid, and this remains 
unchanged by variations in gas pressure. Why then, would reductions in ambient pressure have any influence 
on the contact line’s dynamics? 

An interpretation of the aforementioned phenomena, suggested in a number of different contexts, is that 
the reduction in ambient gas pressure, which leads to a corresponding decrease in gas density, increases the 


mean free path in this phase so that non-equilibrium gas dynamics (Chapman & Cowling 1970), also referred 
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Figure 1: Sketch illustrating the influence of gas flow on the drop impact phenomenon. The behaviour of the 
gas, in the regions indicated by block arrows, affects both (a) the pre-impact and (b) the post-impact spreading 
dynamics of the drop. 


to as rarefied gas dynamics (Cercignani, 2000), become important (Marchand et al.. 2012 Benkreira & Ikin 


2010 Duchemin & Josserand, 2012). Specifically, boundary slip is generated in the gas whose magnitude is 
proportional to the mean free path and thus whose influence is increased by reductions in gas pressure. The 
increased slip allows gas trapped in thin films between the solid and the free-surface to escape more easily, so 
that the lubrication forces generated by the gas are reduced and its influence is negated. 

For drop impact phenomena, whose literature has been previously reviewed ( Rein[ 1993[ Yarin, 2006), 
alterations in the lubrication forces could be important either (a) in the thin air cushion that exists under the 


drop before it impacts the solid (Bouwhuis et al. 2012 Thoroddsen et al.. 2005) and/or (b) in the gas film 


formed in front of the advancing contact line after the drop has impacted the solid (Riboux & Gordillo 
see Figure [I| Theoretical work has also focussed on either (a) the pre-impact stage ([Smith et ^ 2003 


(Schroll et al. 


Shikhmurzaev 


2014), 


Mani 


et at. 2010 Mandre & Brenner 2012) or (b) the post-impact spreading phase where, with a few exceptions 


2010), the gas is usually considered passive (Bussmann et aZ. 2000 Eggers et at. 2010 Sprittles & 


2012a I. Distinguishing how each of these processes influences the spreading dynamics of the drop 


is a complex problem, both from a theoretical and experimental perspective, due to the inherently multiscale 
nature of the phenomenon. As a result, establishing how reductions in ambient pressure suppress splashing 


remains an area of intensive debate (Driscoll & Nagel 2011 Kolinski et al. 2012 de Ruiter et al. 2012) with 


recent experimental results on the impact stage particularly interesting (de Ruiter et al. 2015 Kolinski et al. 


2014). 


In contrast to drop-based dynamics, for the dip-coating flow investigated in Benkreira & Khan (2008) 


mechanism (a) does not exist, as up until the point of wetting failure, the liquid remains in contact with the 
solid and the motion is steady. Therefore, studying this process allows us to isolate the influence of changes in 
ambient pressure on the behaviour of a moving contact line, i.e. effect (b), without any of the aforementioned 
complications that are inherent to the drop impact phenomenon. This is the path that will be pursued in this 
work with drop dynamics only re-considered again in ^ where, having characterised fully (b), our results enable 
us to infer some conclusions about drop splashing phenomena. 


1.1 Experimental Observations in Dip Coating 


The experimental setup constructed in Benkreira &: Khan] ( |2008[ ) is for dip-coating, where a plate/tape is con¬ 
tinuously run through a liquid bath (Figure and gradually increased in speed until air is entrained into the 
liquid. The speed at which this first occurs is called the ‘maximum speed of wetting’ and will be denoted by Uf, 


where stars will be used throughout to denote dimensional quantities. In contrast to previous works (Burley 


& Kennedy 


1976 


Blake &: Shikhmurzaev 


2002), the ambient gas pressure P* can be reduced from its atmo- 


spheric value by up to a factor of a hundred by rnounting this apparatus inside a vacuum chamber. Silicone oils 
of differing viscosities are used as the coating liquid, because their relatively low volatility (compared to, say, 
water-glycerol mixtures) ensures no significant evaporation when the ambient pressure is reduced. This work 
is extended in Benkreira & Ikin (2010) to consider also a range of different ambient gases (air, carbon dioxide 
and helium). The key results from Benkreira & Khan (2008) can be summarised as follows: 


• For a given liquid-solid-gas system, Uf can be increased by reducing P*. 
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Figure 2: Sketch of the dip-coating flow configuration showing a free-surface formed between a liquid (below) 
and gas (above) which meets a solid surface at a contact line, which is a distance Ay below the height of the 
surface (y = 0). The inset shows that the free-surface meets the solid at a contact angle 9d but that a finite 
distance away has an apparent angle 9m at an inflection point in the surface profile where the film has height 


H. 


• There is very little variation in U* until P* is approximately a tenth of its atmospheric value. Further 
reductions in P* lead to sharp increases in U*. 

• For sufficiently low P*, it is possible for high viscosity liquids to have a larger value of U* than lower 
viscosity ones (the opposite of what occurs at atmospheric pressure). 


In addition to these findings, in Benkreira & Ikin (20101 it is shown that: 


• For a given value of P*, using a gas with a larger mean free path increases U*. 


In Benkreira & Ikin (2010), the aforementioned experimental observations are assumed to be caused by 
non-equilibrium effects in the gas which become significant at sufficiently low pressures. These effects, most 
notably slip at the gas film’s boundaries, are interpreted in terms of an ‘effective gas viscosity’ whose magnitude 
is reduced with decreases in ambient gas pressure. Using appropriate values for this parameter, obtained from 


measurements of the resistance of a gas trapped between parallel oscillating plates (Andrews & Harris, 19951, 


allows the data to be collapsed onto a single curve (Figure 9b in Benkreira & Ikin (20101) and highlights the 
influence of non-equilibrium effects. 

Although the effective viscosity may be a useful concept when trying to qualitatively interpret the experi¬ 
mental findings, as an integral property of the system it will not enable us to formulate a local fluid mechanical 
model for this problem. In particular, the dependence of the effective gas viscosity on ambient gas pressure 
which is used in Benkreira & Ikin (2010) was obtained in Andrews & Harris (1995) for a gas film of constant 
height 2/im and will be different for other film profiles. In reality, the gas film’s height will vary all the way 
from the contact line to the far field, where a film can no longer be defined so that the notion of an effective 
viscosity no longer makes sense. Instead, we must look to capture the local physical mechanisms which alter 
the gas film’s dynamics at reduced pressures. 


1.2 Mathematical Modelling of Dynamic Wetting Phenomena 

Consider the various levels of description that can be used to describe dynamic wetting phenomena depending 
on how the gas phase’s dynamics are accounted for, with the simplest case first. 

1.2.1 Dynamic Wetting in a Passive Gas 

In many cases, the viscosity and density ratios indicate that the gas can be treated as dynamically passive. 
Then, mathematical models for the liquid’s dynamics must address two key aspects. Firstly, the classical fluid 
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mechanical model, with no-slip at fluid-solid boundaries, has no solution (Huh & Scriven, 1971 Shikhmurzaev 
20061. This is the so-called ‘moving contact line problem’. And, secondly, the contact angle, which is seen 


to vary at experimental resolutions (Hoffman 1975 Blake & Shikhmurzaev 20021, must be prescribed as a 


boundary condition on the free-surface shape. A review of the various classes of models proposed to address 


these issues can be found in Ch. 3 of Shikhmurzaev (2007). 


To overcome the moving contact line problem, initially identified in Huh & Scriven (19711, it is often assumed 


that some degree of slip occurs at the fluid-solid boundary (Dussan V 1976). Methods for predicting the ‘slip 


length’ of an arbitrary liquid-solid interface are not well developed but estimates suggest it is in the range 


1-10 nm for simple liquids (Lauga et al. 2007). 


The treatment of the dynamic contact angle remains a highly contentious issue. The main question is 


whether or not the experimentally observed variation in the ‘apparent angle’ (Wilson et al.. 2006) is caused 


entirely by the bending of the free surface below the resolution of the experiment, i.e. by the ‘viscous bending’ 
mechanism quantified in jCox ( |1986 ) and measured in Dussan V et a7| ( [l991 ); Rame & Garoff (1996), or whether 
actual contact angle (sometimes referred to as the ‘microscopic angle’) also varies. In either case, the apparent 
angle depends on the capillary number, which is the dimensionless contact line speed, but in the former case 


the actual angle remains fixed whilst in the latter this angle varies, as in the molecular kinetic theory (Blake & 


Haynes 1969). In more complex models, the angle can also depend on the flow field itself, as in the interface 


formation model (Shikhmurzaev, 2007). 


Models in which the contact angle is considered constant have been popular in recent years for investigating 


air entrainment phenomena and have produced a number of promising results (Vandre et al.. 2012 2013 2014 


Marchand et al. 

2012 

Chan et al. 

2013 


surface becomes significant as the point of air entrainment is approached and cannot be neglected. What 
remains unclear is whether this mechanism is sufficient to account for the experimentally observed dynamics 
of the apparent contact angle, particularly given that precise measurements of the magnitude of slip on the 
liquid-solid interface are lacking. These difficulties have stimulated fundamental research into the contact line 


region’s dynamics using molecular dynamics techniques (De Coninck & Blake 2008) but, as can be seen from a 
recent collection of papers (Velarde 2011) and review articles (Blake 2006 Snoeijer fc Andreotti[ 2013), intense 
debate remains about the topic. 

The focus of this work is to characterise the effects which the gas dynamics can have on dynamic wetting 
phenomena and, consequently, the simplest possible model for the other aspects will be used. Specifically, 
slip at the liquid-solid interface will be captured using the Navier condition (Navier 1823), as commonly used 


in dynamic wetting models (Hocking, 1976 Huh & Mason, 1977), and the contact angle will be taken as a 


parameter that is fixed at its equilibrium value. The result is a ‘conventional model’ and a comprehensive 


review of all such models which can be ‘built’ is given in Ch. 9 of Shikhmurzaev (1997). More complex models 
can easily be constructed on top of the model used here, but in the present context would distract from the 
main focus of this work. 


1.2.2 In a Viscous Gas 

It is now well recognised that at sufficiently high coating speeds a thin film of gas forms in front of the moving 


contact line, as observed experimentally by laser-Doppler velocimetry in Mues et al. (1989). The thin nature 


of the film means that despite the large viscosity ratio, lubrication forces enhance the influence of the gas’ 
dynamics so that they cannot be neglected. As coating speeds are increased, the resistance this film creates 


to contact line motion eventually becomes sufficiently large that it is entrained into the liquid (Vandre et al. 


2013 Marchand et al ., 2012). Therefore, to study the phenomenon of air entrainment and wetting failures, any 


theory developed must account for the dynamics of the viscous gas as well as the liquid. 

In formulating a model for the gas’ dynamics one again encounters the moving contact line problem, namely 
that if no-slip is used then the contact line can’t move. Consequently, slip has also often been applied on the 
solid-gas boundary with a magnitude, measured by the slip length, either explicitly stated to be the same as in 
the liquid (Vandre et al. 2013) or implicitly assumed to be so (Cox 1986). In other words, slip is usually used 


as a mathematical fix to circumvent the moving contact line problem, rather than as a physical effect whose 
magnitude can significantly alter the flow. 


1.2.3 In a Non-Equilibrium Gas 


Experimental observations in Marchand et al. (2012), where the film entrained by a plunging plate was measured 


using optical methods, suggest that for relatively viscous liquids the characteristic film thickness is rather small, 
in the range H* ^ 1-15 /rm. For flows where there is an external load on the film, such as in curtain coating, 
the film’s thickness is likely to be even smaller. Consequently, the film’s dimensions become comparable to the 
mean free path in the gas £* so that the Knudsen number Kuh = £*/H*, characterising the importance of 
non-equilibrium gas effects, becomes non-negligible, particularly at reduced pressures. For example, a factor of 
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Figure 3: The effect of Maxwell-slip on the gas flow. When A = 0 slip is present at the gas-solid boundary but 
turned off at the gas-liquid boundary, so that Ug = Ug, whilst for A = 1 it is present on both interfaces. 


ten decrease in the ambient pressure of air, achieved in both Benkreira & Khan (2008) and Xu et al. (20051, 
gives a value of = 0.7 /rm so that Kuh ~ 0.1. In this case, one may expect that non-equilibrium effects in 
the gas film will have a significant impact on its dynamics, and we must consider how this can be incorporated 
into our dynamic wetting model. 

Non-equilibrium effects first manifest themselves at the gas’ boundaries and can be accounted for by relaxing 
the no-slip condition to allow for boundary slip whilst continuing to use the the Navier-Stokes equations in the 
bulk. This approach remains accurate for Kuh < 0.1, during which gas flow is said to be in the ‘slip regime’, 
after which Knudsen layers will begin to occupy the entire gas film, so that the scale separation between 
boundary and bulk effects no longer exists. 


In Maxwell (1879), a slip condition was derived for the ‘micro-slip’ at the actual boundary of the gas whose 
mathematical form remains the same, when Knu ^ 1, if ‘macro-slip’ across a Knudsen layer of size ^ £* is 


also accounted for, see §3.5 of Cercignani (2000). We shall henceforth refer to the discontinuity in tangential 


velocity as ‘Maxwell slip’ which at an impermeable boundary is generated by the stress acting on that interface 
from the gas phase. For two-dimensional isothermal flow, the component of velocity m*, where subscripts g will 
denote properties of the gas, adjacent to a stationary planar surface located at y = 0 is, in dimensional terms, 
given by: 

du* 


at 


dy* 


= Un 


( 1 ) 


where the coefficient of proportionality a can depend on both the properties of the micro-slip, which will 
depend on the fraction of molecules which reflect diffusively and those which a purely specular, and the macro¬ 
slip generated across a Knudsen layer. In practise, unless surfaces are molecularly smooth we will have a ^ I 
(Agarwal & Prabhu 2008 Millikanj |1923 Allen & Raabe, 1982) and, for simplicity, we will henceforth take 
a = 1. 

The Maxwell slip condition ([^ has the same form as the slip conditions which have often been used in 
dynamic wetting flows to circumvent the moving contact line problem in the gas but occurs naturally here from 
the consideration of non-equilibrium gas effects. Notably, this means that accurate values of the slip length 
can easily be obtained by knowing the mean free path in the gas, in contrast to the slip length in the liqnid 
which is notoriously difficult to determine. Furthermore, the approach taken here suggests that slip will also 
be present at the gas-liquid boundary (i.e. the free-surface), as demonstrated by experiments dating back to 
Millikan (1923). This effect is shown in Figure]^ where the parameter A is used to change between the usual 
condition of continuity of velocity across the gas-liquid interface (A = 0) and Maxwell-slip (A = 1). 

Slip at the gas-liquid boundary has been accounted for in some models for the gas film in drop impact 


(Duchemin & Josserand 2012 Mani et al.. 2010), and in great detail for drop collisions (Sundararajakumar & 


Koch 1996), but is yet to be routinely incorporated into dynamic wetting models. In most cases where slip 


on the gas-solid interface is accounted for, the same effect on the gas-liquid surface is not considered (Vandre 


et al. 


et al. 


2013 


Riboux & Gordillo, 2014). The only progress in this direction has been achieved in Marchand 


(2012), where non-equilibrium effects are accounted for by modifying the viscosity in the gas. However, 


this is achieved by using an empirical expression for the effective viscosity which suffers from the same issues 
mentioned previously, namely a reliance on the lubrication approximation that does not enable one to fully 
formulate a uniformly-valid local fluid mechanical problem. 
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In order to apply Maxwell-slip at curved moving surfaces, the tensorial expression is required (Lockerby 
et al. 2004). Specifically, for gas flow adjacent to a surface which moves with velocity u* and has normal n 


pointing into the gas, we have 


t n ■ (Vu* -f [Vu*]^) . (I 


nn) =u*|| -ufj. 


( 2 ) 


where the subscripts || denote the components of a vector tangential to the surface which can be obtained for 
an arbitrary vector a using a|| = a • (I — nn), where I is the metric tensor of the coordinate system. 

If the gas is modelled as being composed of ‘hard spheres’, then the mean free path is related to the 
temperature T* and ambient pressure P* by P = k*sTy{V27TP*iP^J^) where is the molecular diameter 
and kg is Boltzmann’s constant. Notably, if local variations in gas pressure about P* are small, which will be 
the case here, then for a given ambient pressure i* is fixed throughout the entire gas phase. In practise, when 
considering the dependence of the mean free path of a gas on the ambient pressure at a fixed temperature it is 
convenient to use the expression 



P* 

g,atm 


9 


1/^9 


(3) 


where Pg is the gas pressure normalised by its atmospheric value. At atmospheric pressure, air has a mean free 
path P’a.tm — other commonly encountered gases take similar values ~ 0.1 /rm. 

Although the Knudsen number in the gas film Kuh is of most interest, this can only be obtained once our 
computations have allowed us to find H*. Before then, the Knudsen number based on our characteristic scales 
Kn = P/L*, where L* ~ 1 mm is the capillary length, must be used. Noting that the dimensionless mean free 
path at atmospheric pressure is £atm = i f^e resulting expression for Kn takes the form 


Kn=^-^ (4) 

^9 

SO that the explicit dependence on pressure reduction becomes clear. 

1.3 Overview 

This work will answer the following questions: 


• Can increases in Maxwell slip, caused by experimentally-attainable reductions in ambient gas pressure, 
significantly raise the maximum speed of wetting? 

• Does slip on the liquid-gas interface, usually neglected in this class of flows, have a substantial effect on 
the dynamic wetting process? 


Will the magnitude of the Knudsen number based on the gas film’s height Knn be small enough for the 
gas flow to remain in the ‘slip regime’? 


Is it possible to collapse the data from Xu et al. (2005) onto a master curve by assuming that non¬ 
equilibrium gas effects are responsible for splash suppression in drop impact? 


The paper is organised as follows. In ^ the dip-coating problem will be formulated and dimensionless 
parameters will be identified. The mathematical model for this process will be seen to be inherently multiscale, 
with flow dynamics on the nano-, micro- and mm-scale caused by, respectively, slip in the liquid, slip in the gas 
and capillarity. The computation of such multiscale flows is challenging and neglecting to resolve all scales in the 


problem yields inaccurate results (Sprittles & Shikhmurzaev 

2012c 

). In ^ the application of a computational 

framework developed in 

Sprittles & Shikhmurzaev 

(2012cl 120131) to this problem is described, which resolves all 


scales, and our method for defining a maximum speed of wetting is shown. Benchmark computations performed 
at atmospheric pressure are provided in the Appendix. In ^the effects of the non-equilibrium gas dynamics on 
a base state are analysed and the gas flow is described in a lubrication setting. After, in ^ a parametric study 
of the new system is conducted which allows us to identify the importance of the various competing physical 
mechanisms in the experimentally-verifiable regimes before, in ^ directly comparing our computations to 
experimental data. The implications of these findings for the drop impact phenomenon are considered in 0 
Finally, in ^ potential avenues of enquiry for future theoretical and experimental development in this area are 
highlighted. 
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2 Problem Formulation 


Consider the flow generated by the steady motion of a smooth chemically-homogeneous solid surface which is 
driven through a liquid-gas free surface in a direction aligned with gravity at a constant spe ed U* (see Figure[^ . 

A characteristic length scale for this problem is given by the capillary length L* = yjcr*/{{p* — Pg)g*) — 

yj<7*/{p*g*)^ where cr* is the surface tension of the liquid-gas surface, which is assumed constant; g* is the 
acceleration due to gravity; and p*, Pg are the densities in the liquid and gas, respectively, with p'^jp* 1 in 
the liquid-gas systems of interest so that the capillary length can be determined solely from the liquid’s density. 
A scale for pressure is p*U*/L*, where p* is the viscosity of the liquid and p* will be its value in the gas. 


2.1 Bulk equations 

For compressibility effects in the liquid and gas to be negligible we require that the Mach number Ma = U* jc* 
in each phase remains small, where c* is the speed of sound. Given that U* <1 ms“^ and c* > 100 ms“^ in all 
cases considered, we have Ma < 0.01. However, this is not a sufficient condition for the flow to be considered 
incompressible, especially when there are thin films of gas present along which pressure can vary significantly, 
see §4.5 of Gad-el-Hak (2006). Here, we will assume that the flow in both phases can be described by the 
steady incompressible Navier-Stokes equations and in §4.5| the validity of this assumption will be confirmed 
from a-posteriori calculations. Therefore, we have 


V • u = 0, Ca u • Vu = Oh^ (V • P -I- g/Ca), 

p Ca Ug • Vug = Oh^ {V - Pg + p g/Ca). 
where the stress tensors in the liquid and gas are, respectively. 


V • Ug = 0, 


(5) 

( 6 ) 


P = —pi + Vu-l-(Vu) 


and 


P9 = -PsI + A^ Vug-l-(Vug) 


(7) 


Here, u and Ug are the velocities in the liquid and gas; p and Pg are the local pressures in the liquid and gas, 
in contrast to uppercase P’s which are used for the ambient pressure in the gas; Ca = p*U*/a* is the capillary 
number based on the viscosity of the liquid; the viscosity ratio is p = p*/p*', the ratio of gas to liquid densities is 
p = p*g/p*] and g is a unit vector aligned with the gravitational field. The Ohnesorge number Oh = p*/yjp*a*L* 
has been chosen as a dimensionless parameter, instead of the Reynolds number Re = Ca/Oh^ = p*U*L*/p*, as 
for a given liquid-solid-gas combination Oh will remain constant, i.e. it will depend solely on material parameters 
and be independent of contact-line speed. 


2.2 Boundary Conditions at the Gas-Solid Surface 

Conditions of impermeability and Maxwell slip give 

Ug • ris = 0, Kn • Pg • (I - ngn^) = Ug|| - Uy, (8) 

where U is the velocity of the solid and the normal to the solid surface is denoted as . 


2.3 Boundary Conditions at the Liquid-Gas Pree-Surface 

On the free-surface, whose location must be obtained as part of the solution, for a steady flow the kinematic 
equation and continuity of the normal component of velocity give that 


u • n = Ug • n = 0, (9) 

where n is the normal to the surface pointing into the liquid phase. These equations are combined with the 
standard balance of stresses with capillarity that give 


n • (P - Pg) • (I - nn) = 0, n • (P - Pg) • n = V • n/Ca. 


( 10 ) 


Equations ([^ and (10) are usually combined with a condition stating that the velocity tangential to the 
interface is continuous across it Ug|| = uy . However, when Maxwell slip is accounted for this is replaced with 


— A Kn n • Pg • (I — nn) = Ugy — uy, 


( 11 ) 


where the minus sign occurs because the normal points into the liquid. By setting A = 0, Maxwell slip can be 
‘turned off’ at the gas-liquid interface, see Figure]^ as has been the case in previous dynamic wetting works 
which only consider slip at the solid boundary. 
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2.4 Boundary Conditions at the Liqnid-Solid Snrface 

The standard conditions of impermeability and Navier slip give 

u-118 = 0, 4 n • P • (I - rign^) = U|| - U||, (12) 

where Ig = is the dimensionless slip length based on the (dimensional) slip length at the liquid-solid 

interface I*. This parameter has no relation to the slip coefficient in the gas which was based on the mean free 
path £*. In other words, although Maxwell slip and Navier slip have the same mathematical form, their physical 
origins differ and thus there is no reason to expect their coefficients to have similar magnitudes. 


2.5 Liquid-Solid-Gas Contact Line 

Equation ( |Io| ) requires a boundary condition at the contact line which is given by prescribing the contact angle 
9d satisfying 

n • rig = — cos0d. (13) 


Here, our focus is on understanding the effects of the gas dynamics on the dynamic wetting flow so that we 
will take the simplest option of assuming 64 to be a constant. More complex implementations, such as those in 


Sprittles & Shikhmurzaev (2013), can be considered in future research. 


2.6 ‘Far Field’ Conditions 

Dip coating is usually conducted in a tank of finite size whose dimensions are large enough not to alter the 


dynamics close to the contact line, as confirmed in Benkreira & Khan (20081. The influence of system geometry. 


or ‘confinement’, on the dynamic wetting process are well known (Ngan & Dussan V 1982) and, in particular. 


have been used to delay air entrainment (Vandre et al. 2012), but these effects are not the focus of this work. 


Therefore, the ‘far field’ boundaries of our domain, which are assumed to be no-slip solids u = Ug = 0, are set 
sufficiently far from the contact line that they have no influence on the dynamic wetting process. Where the 
free-surface meets the far-held it is assumed to be perpendicular to the gravitational held, i.e. ‘hat’, so that 
n • g = 1. In practise, setting the far held (dimensionally) to be twenty capillary lengths from the contact line 
in both the x and y directions was found to be sufficiently remote. 


3 Computational Framework 


Dynamic wetting and dewetting hows have often been investigated using the simplihcations afforded by the 


lubrication approximation (Voinov 1976 Eggers 2004). However, for the two phase how considered here, this 


approximation cannot simultaneously be valid in each phase, i.e. the liquid and gas phases can’t both be thin 
hlms. As shown in recent papers (Marchand et al. 2012 Vandre et al.. 2013), plausible extensions in the spirit 


of the lubrication approach can be attempted, but their accuracy can only be ascertained from comparisons 
to full computations and is unlikely to be acceptable at larger capillary numbers where the contact line region 
cannot be considered as ‘localised’. 

At present, to obtain accurate solutions to the mathematical model formulated in one is left with no 
option other than to use computational methods. 

Computationally, resolving all scales in the problem is particularly important to ensure that our contact line 


equation (13) is applied to the ‘actual’ or ‘microscopic’ angle 9^, rather than the ‘apparent one’ 9, 
i.e. the angle of the free surface at a finite distance s from the contact line. Although the apparent angle can. 


app ^app (^)5 


in certain circumstances, be related to the actual angle (Cox 1986), this will not be true at the high capillary 


numbers observed in coating flows, particularly when non-standard gas dynamics are included. 

Here, gas dynamics will be built into a finite element framework developed in [Sprittles fc Shikhmurzaev] 
(2012c 2013) to capture dynamic wetting problems and then extended to consider two-phase flow problems in 
our work on the coalescence of liquid drops (Sprittles & Shikhmurzaev, 20146). Notably, this code captures 


all physically-relevant length scales in the problem, from the slip length I* ~ 10 nm up to the capillary length 
L* ^ 1 mm, meaning that at least six orders of magnitude in spatial scale are resolved. As a step-by-step 
user-friendly guide to the implementation has been provided in Sprittles & Shikhmurzaev (2012c), only the 


main details will be recapitulated here alongside some aspects which are specific to the current work. 


The code uses the arbitrary Lagrangian Eulerian scheme, based on the method of spines (Ruschak 1980 


Kistler & Scriven 1983), to capture the evolution of the free surface in two-dimensional or three-dimensional 
axisymmetric flows. Both coating flows, which are often time-independent, as well as unsteady flows such 


as drop impact (Sprittles & Shikhmurzaev 20126) and drop coalescence (Sprittles & Shikhmurzaev 2012a 


20146|a ) have been considered. It has been confirmed that the code accurately captures viscous, inertial and 


capillarity effects, even when the mesh undergoes large deformations which inevitably occur when Ca ~ 1. 
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The mesh is graded so that small elements can be used near the moving contact line to resolve the slip length 
whilst larger elements are used where only scales associated with the bulk flow are present. Consequently, the 
computational cost is relatively low so that even for the highest resolution meshes used in this work, for a given 
liquid-solid-gas combination the entire map of, say, Ca vs contact line elevation (Ay), can be mapped in just a 
few hours on a standard laptop. 

Computations at high Ca, roughly those for which Ca > 0.1, are notoriously difficult as (a) the free-surface 
becomes increasingly sensitive to the global flow and (b) regions of high curvature at the contact-line require 
spatial resolution on scales at or below the slip length Ig- Both factors hinder the ease at which converged 
solutions can be obtained and these issues are compounded when Oh is small, indicating the importance of 
non-linear inertial effects. For all parameter sets computed, converged solutions can always be obtained for 
Ca < 2, which compares very favourably to previous works, so that this is chosen as an upper cut-off for the 
results of our parametric study. 


3.1 Benchmark Calculations for the Maximum Speed of Wetting 


In coating flows, one would like to know for a given liquid-solid-gas combination the wetting speed U* at 
which the contact line motion becomes unstable so that gas is entrained into the liquid either through bubbles 


forming at the cusps of a ‘sawtooth’ wetting-line (Blake & Ruschak 19791 or by an entire film of gas being 


dragged into the liquid (Marchand et al. 2012). This is referred to as the ‘maximum speed of wetting’ and, in 
dimensionless terms, is represented by a critical capillary number Cac = ^*U*/a*. Our method for calculating 
Coc is explained in the Appendix alongside benchmark calculations for its value which are compared to the 
results of! 
results in 


Vandre et al.\ (|2013[) across a range of viscosity ratios. Excellent agreement is obtained between the 


Vandre et al. (2013) and those found here. As these calculations distract from the main emphasis of 


this work, but could be useful as a benchmark for future investigations in the field, they are provided in the 
Appendix. 


4 Effect of Gas Pressure on the Maximum Speed of Wetting: Anal¬ 
ysis of a Base State 

Having demonstrated the accuracy of our code and shown how the critical capillary number Cac is calculated, 
we now investigate the relation between the gas pressure Pg and Ca^ Here, we will analyse a base state, varying 
only gas pressure and A, in an attempt to deepen our understanding of the influence of non-equilibrium gas 
effects on the dynamic wetting process before, in performing a parametric study of the system of interest. 

To ensure we are studying an experimentally attainable regime, we use the system considered in |Benkreira| 
& Khan (2008) to provide a base state about which our parameters can be varied. In this work, a range of 
silicone oils of different viscosities /i* = 20-200 mPa s are used whilst the density p* = 950 kg m“^, surface 
tension a* = 20 mN m“^ and equilibrium contact angle 9c = 10° remain approximately constant. Then the 


characteristic (capillary) length L* = 1.5 mm and Oh = 6 x 10 (mPa s) ^ depends only on the viscosity 

p*. 

A base case is chosen by taking p* = 50 mPa s so that Oh^ = 0.3, where base state values will be denoted 
with a subscript 0. Taking air as the displaced fluid with p* = 18 /iPa s, which remains independent of the 
ambient pressure (Maxwell 1867), gives a viscosity ratio of po = 3.6 x 10“"^. The density of the gas p* depends 
on the ambient gas pressure, with its maximum value of p* = 1.2 kg m“^ at atmospheric pressure making 
Pmax = 1-3 X 10“^. However, for all calculations performed, the value of p has a negligible effect on Cac when 
P ^ Pmax and so henceforth, without loss of generality, we take p = 0. 

The slip length of the liquid-solid interface is fixed at I* = 10 nm, which is well within the range of 
experimentally observed values (Lauga et al. 2007), so that the dimensionless parameter C = 6.8 x 10“®. 


Assuming the gas is air, we have Kn = 4.8 x 10 ^/Pg so that the only parameter which remains to be specified 
is A characterising whether or not there is slip at the gas-liquid boundary. 


4.1 Effect of Gas Pressure 

In Figure]^ the effect on Cac of reducing Pg is computed for A = 0 and A = 1. In each case, for Pg > 0.1 there 
is only a slight increase in Cac from its base value of 0.47, which is independent of A. For Pg < 0.1, dramatic 
changes in Cac are observed, whose form depends on A, so that Pg « 0.1 appears to be the point at which 
non-equilibrium effects become important. 

Notably, as can be seen most clearly in the inset of Figure]^ for A = 0 as Pg 0 we have Cac —> 0.64 
whereas for A = 1 computations suggest that Cac —> oo. In fact, in the latter case, once a critical pressure Pg^c 
is approached, Cac appears to increase without bound. For A = 1 this occurs at Pg^c = 9 x 10“^; just above 
this value, at Pg c = 10“^, we have Cac = 1-09 whilst for all Pg < A,c we find Cac > 2. 
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Figure 4: The dependence of the critical capillary number Cuc on the gas pressure Pg, for our base parameters, 
with curve 1: A = 0 (in black) and 2: A = 1 (in blue). The inset shows how curve 1 (A = 0) tends to a finite 
value of Coc as Pg —>■ 0 whilst curve 2 appear to increase without bound. 


It is impossible to show rigorously that Coc —>■ oo as Pg —>■ Pg c without either the computation of higher 
values of Cqc, which would allow some scaling to be inferred, or, ideally, the development of an analytic 
framework for this problem. Unfortunately, in the former case (see §3) computational constraints permit us 
from reaching Cac > 2 whilst in the latter we are confined by the fact that all analytic development, in particular 
those of Cox (1986), are only valid for small Cac- Despite this, the rapid increase of Cuc as Pg,c is approached 


is striking and Coc = 2 is still a large value in the context of dip-coating phenomena. 

To summarise, it has been shown that when Maxwell slip is accounted for on both the gas-solid and gas- 
liquid interfaces the maximum speed of wetting appears to be unbounded as the ambient pressure is reduced, 
whereas if slip on the latter surface is neglected (^ = 0), the maximum speed remains finite. In other words, 
the error associated with neglecting Maxwell slip on the free-surface is extremely large at reduced pressures and 
appears to be infinite in the limit Pg —>■ 0. 


4.2 Flow Kinematics 

In Figure]^ streamlines of the flow near the contact line are shown at Ca = 0.4 in three different regimes: (a) at 
a pressure where non-equilibrium effects are weak Pg = 0.14 with 4. = 0, (b) Pg = 0.14 with A = 1 and (c) at a 
pressure where non-equilibrium effects are becoming influential Pg = 0.014 with A = 1. Given that L* ~ 1mm 
for typical liquids, dimensionally the scale in Figurej^is a few microns. In all cases, the flow of the liquid, which 
is below the free-surface (represented by a thick black line), remains virtually unchanged, with the motion of 
the solid, located at x = 0, dragging liquid downwards which, to conserve mass, is continually replenished from 
above. If one notes that the free-surface meets the solid at an equilibrium contact angle of 10°, and yet any 
apparent angle defined on the scale seen in Figure would be obtuse, it is clear that at this relatively high 
capillary number there is significant deformation of the free-surface on a scale below what is visible here. 

In all three cases, the flow of liquid, which seems immune to the gas’ dynamics, results in an almost identical 
velocity tangential to the free-surface on the liquid-facing side of this interface, as shown from curves la, 2a, 3a 
in Figure]^, where the velocities tangential to the gas-solid and liquid-gas interfaces have been plotted. 

At the top y = 2.5 of the figures, the gas flow field is qualitatively similar in all cases. The motion of both 
the liquid and the solid drives gas towards the contact line region which, to ensure continuity of mass, results 


in a ‘split ejection’ type flow, as observed in immiscible liquid-liquid systems (Dussan V & Davis, 1974 Dussan 
|V 1977), with an upwards flux of gas through the middle of this domain. For A = 1 the split ejection flow is 
maintained right up to the contact line (Figure [^,c); however, for A = 0 a flow reversal occurs on the solid-gas 
interface at y « 1 so that the direction of the gas flow actually opposes the solid’s for y < 1 (Figure §1). 
The flow reversal can clearly be seen from curve 1 in Figure for s < 10 with a substantial minimum of 
U 2 t = —0.3, and can also be seen in previous works, such as Figure 10b of Vandre et al. (2013). 


The asymmetric gas flow observed for A = 0 appears because in this case the gas’ velocity is continuous 
across the liquid-gas free-surface, see curves la and lb in Figure [§i, but is allowed to slip across the solid 
boundary according to Maxwell’s condition, see curve 1 in Figure [^. This results in the liquid driving a gas 
flow towards the contact line which, when combined with the requirement of mass conservation, leads to such 
a large build of tangential stress at the solid, that the slip across the interface, determined by pL is sufficient 
to reverse the flow direction. This effect will be further considered, in a lubrication setting, in §4.4| 
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Figure 5: Streamlines computed at Ca = 0.4 for different ambient gas pressures Pg both with {A = 1) and 
without {A = 0) Maxwell slip at the free-surface. In (a): {Pg,A) = (0.14,0), (b): (0.14,1) and (c): (0.014,1). The 
local ‘zoomed in’ coordinates used are x = x x 10^ and y = {y — yc) x 10^, where yc is the contact line position, 
and streamlines emanate from equally spaced points across y = 2.5 with the dividing streamline dashed. 
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Figure 6: Tangential velocities as a function of distance from the contact line s along the (a) liquid-gas and (b) 
gas-solid boundaries computed at Ca = 0.4 with curve 1: {Pg,A) = (0.14,0), 2: (0.14,1) and 3: (0.014,1). The 
velocity tangential to the free-surface and pointing towards the contact line is uu, with subscripts ‘a’ for the 
liquid-facing side of the interface and subscripts ‘b’ for the gas side. The tangential velocity on the gas-solid 
interface pointing towards the contact line is U 2 t and the speed of the solid is 1. 
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Figure 7: Free-surface shapes obtained in our base case, with Pg = 0.014, with curve 1: Ca = 0.4 (in black), 
2: Ca = 0.6 (in red) and 3: Ca = 0.8 (in blue). The inset shows how the height of the gas film h varies as a 
function oi y — yc showing that for 10“"^ < y — Vc < 10“^ it can be considered in a lubrication setting. 


When A = 1, slip is also allowed at the gas-liquid interface, see curves 2a and 2b in Figure]^, so that, as 
can be seen from curve 2 in Figure [^, there is no flow reversal on either of the gas’ boundaries. This more 
symmetric flow field remains when the gas pressure is decreased to Pg = 0.014, see Figure]^ and curves labelled 
3 in Figurej^ Although the flow remains qualitatively the same, one can clearly see from Figurej^that lowering 
the pressure significantly increases slip at the gas’ boundaries. For example, at s = 10“^, for Pg = 0.14 the 
velocity on the solid U 2 t = 0.65 whereas for Pg = 0.014 the velocity is just U 2 t = 0.16. Therefore, less gas is 
dragged into the contact line region when (a) slip is accounted for on the gas-liquid boundary [A = 1) and/or 
(b) slip is increased through reductions in pressure. Both of these mechanisms reduce the gas’ resistance to 


contact line motion, as shown analytically in 14.4 and contribute to postponing the point of wetting failure. 


4.3 Characteristics of the Gas Film 

At high capillary numbers there is a region between the contact line and the far field in which the gas domain 
is a ‘thin film’. This region is sufficiently far from the contact line region where the free-surface bends from its 
equilibrium angle of 10° towards an obtuse apparent angle and near enough to the contact line that gravitational 
forces have not started to flatten the free-surface. For our base case, for Ca > 0.4, which includes all values 
of Coc, this region exists between 10“"* < y — yc < 10“^ as one can see from FigureSpecifically, the main 
figure clearly shows that the gas film is ‘thin’ on the scale of around y — j/c ~ 0.1 whilst the inset shows on a 
logarithmic plot that this behaviour continues down to y — yc ^ 10“^. 

The dynamics of the gas film are key to the air entrainment phenomenon so that one may expect that 
reductions in ambient pressure will only start to influence Cac once Maxwell slip on the boundaries is large 
enough to affect the gas’ flow characteristics. This happens when the (dimensionless) slip length Kn = 4.8 x 
10“®/Pg becomes comparable to the height of the gas film x = h{y). This results in a. function Knioc{h) = Kn/h 
which characterises the importance of non-equilibrium effects as one goes along the film. How though, should 
we define a particular position along the film h = H on which we can base a local Knudsen number Knu"^. 


Research into dynamic wetting/dewetting phenomena (Dejaguin & Levi 1964 Vandre et al. 2013 Eggers 


20041 has long recognised that wetting failure is strongly linked to the behaviour of the inflection point on 


the free-surface, where the surface’s curvature changes sign; this point can also be used to define an apparent 


contact angle (Tanner 19791. Given that the flow near the inflection point seems to be key to determining Coc, 
this point is a sensible choice for defining the characteristic film height H, see Figure To do so, we must 
extract H from our computations and determine how it depends on Ca and Pg. 

In Figure]^ we can see that although Coc depends on Pg, the value of H is approximately independent of 
Pg, so that H « H{Ca), with smaller Pg simply revealing more of this curve. In fact, a similar curve is obtained 
if the gas flow is ignored altogether {p, — 0), as shown by the dashed line in Figure]^ although this may not be 
the case at higher viscosity ratios. Notably, H decreases rapidly with increasing Ca in agreement with previous 
experiments (Marchand et al. 20121 and simulations (Vandre et al. 2013[ ). Reassuringly, at Ca = Cac the film 
height obtained at atmospheric pressure, H = 4.2 x 10“'^, dimensionally corresponds to a film H* = 1-10 ym 
for typical liquids, which is precisely what has been found experimentally in Marchand et al. ( 2012[ ). 

Enhancements in Coc at reduced pressure can be quantified by the percentage increase in Cac from its 
atmospheric value Cac.atm defined by ACoc = 100(Cac — Cac^atm)/Cac^atm, which is given in Table[^ As well 
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Figure 8: The dependence of the gas film height H, at the inflection point on the free surface, on the capillary 
number Ca. The curve is largely independent of gas pressure Pg, although the lower values of Pg allow more 
of the curve to be obtained. Curves corresponds to 1: Pg = 1 (red), 2: Pg = 10“^ (blue) and 3: Pg = 9 x 10“^ 
(green), with the dashed line obtained for fl = 0. 


ACoc 

Pg 

H 

Kuh 

0 

10% 

20% 

50% 

100% 

> 300% 

1 

0.16 

7 X 10-2 

2.3 X 10-2 
1.2 X 10-2 
< 9 X 10-3 

4.2 X 10-3 

3.1 X 10-3 

2.3 X 10-3 

1.1 X 10-3 

4.1 X 10-“^ 

< 7.1 X 10-3 

0.011 

0.1 

0.3 
1.9 
9.8 
> 75 


Table 1: The pressure reduction Pg required for a given enhancement in capillary number ACuc with the 
corresponding film height H and Knudsen number at the inflection point Knu at this pressure. 


as listing the pressure reduction required, values for H are given from which Kn^ = Kn/PI is calculated. The 
data shows that for Pg ~ 0.1, where significant increases in ACqc begin, corresponds to Knn ^ 0.1. Therefore, 
non-equilibrium effects in the gas, which manifest themselves through Maxwell slip, alter the flow once the slip 
length is around 10% of the gas film’s characteristic height and it is this mechanism which determines when 
variations in Pg can start to affect Cqc- 

Notably, all inflection points calculated in Table [l] fall into the thin film region apart from the entries where 
the critical gas pressure has been passed after which the inflection point gets close to the region of high curvature 
near the contact line. This suggests that the critical pressure Pg_c = 9 x 10“^ is the one which is low enough 
for the inflection point to no longer be located the thin film region. From Table we can see that for Pg < Pg ,. 
we have Knn > 75 so that Maxwell slip is so strong that the gas flow is hardly affected by the motion of its 
boundaries. In this case, the gas effectively offers no resistance to the dynamic wetting process and Cqc is able 
to increase without bound. 


4.3.1 Limitations of the MaLxwell-slip Model 

From Table it is clear that the combined increases in Kn and decreases in PI at reduced pressures both 
contribute to rapid changes in Kuh- Across a reduction in the gas pressure of ~ 10^, Kuh increases from its 
atmospheric value by a factor of ^ 10^. Consequently, most values of Knn in the Table fall well outside the 
‘slip regime’ where non-equilibrium effects can be attributed entirely to the boundary conditions. Therefore, 
attributing non-equilibrium effects entirely to the boundary conditions, via Maxwell-slip, as has been considered 
here is not sufficient to accurately capture its behaviour. This calls into question a number of the results which 
have been obtained with a gas model which goes outside its strict limits of applicability for the smaller values 
of Pg. In particular, is it likely that the behaviour Cqc —>■ oo as Pg —>■ Pg,c, or even as Pg 0, remains robust? 

To address this question we must understand more about the Knudsen layer in regimes where its effects 
can no longer be attributed entirely to the boundary conditions. A useful result for this phenomenon has been 


obtained in Lockerby et al. (20056), for pressure-driven Poiseuille flow in a channel. There, it was shown that 
at a Knudsen number of just Kn = 0.025, the increased mass flow rate above that obtained for no-slip is 
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Figure 9: Vertical velocity v predicted by the lubrication analysis for h = 2.5 x 10 “^ at Pg = 0.14 for the cases 
of A = 0 , curve 1 in red, and A = 1, curve 2 in black. One can see that when A = 0 the lubrication analysis 
predicts the flow reversal at x = 0 whilst this does not occur for A = 1 . These results can be compared to those 
in Figure[^,b by looking at y = 0.3 where h = 2.5 x 10 “"*. 


approximately 70% due to increased slip at the wall whilst 30% is due to non-Newtonian effects in the Knudsen 
layer. This latter phenomenon is not accounted for in our work and suggests that incorporating more complex 
gas dynamics is likely to enhance the effects already observed with Maxwell-slip rather than suppress them. 
Therefore, there is good reason to believe that with the incorporation of more complex gas models the qualitative 
trends observed for Coc will remain whilst quantitatively its value could increase for a given Pg. 

We shall return to these points in more detail in however, in what follows we will continue to use the 
problem formulation outlined in ^ to give some insight into the role that non-equilibrium effects play in the 
gas, despite, in some cases, the model being outside its strict region of applicability. 


4.4 Lubrication Analysis 


Given the importance of the gas film’s behaviour, it is of interest to see if analytic progress in a lubrication 
setting will shed some light on the dynamics of this film. Figurej^ shows that uit, the liquid’s velocity tangential 
to the free-surface is within 10% of uu = 0.85 throughout the thin him region. This is a useful observation, as 
it allows us to make analytic progress by assuming that the liquid approximately provides a constant downward 
velocity V = 0.85 along the free-surface in this region. 

In the lubrication setting, the steady gas how between two impermeable surfaces, at x = 0 and x = h, 
generated by their velocities of magnitude, respectively, v = —1 and v = —V, will, in order to conserve mass 
{Jq V dx = 0), have parabolic form v = a(x^ — h^/3) -|- h{x — /i/2) . The coefficients a,b are obtained by 
applying the boundary conditions at x = 0, h which are, respectively, the Maxwell-slip equations 
Introducing Kriioc = Kn/h, the result is 


and (111. 


-3[l + 2AKnioc + Vi2Knioc + l)] 2(3 - ah^) 

“ fPll + 4^moc{l + A) + 12AKn:fj “ 3h {2Knioc +1)' 

A good test for the lubrication theory is to see whether it is able to predict how the reversal at the solid 
surface observed in Figure]^ for the case of A = 0. To do so requires that v{x = 0) > 0, which for A = 0 can be 
shown to occur when Kriioc > 1/(21^) so that how reversal is indeed possible if Knioc is sufficiently large. For 
the case in Figure [^, using V = 0.85 and Kriioc = 3.4 x 10“^//i, the lubrication analysis predicts how reversal 
for h < hr = 5.8 x 10“"^ with a typical how prohle in this regime shown by curve 1 in Figure]^ The agreement 
with the computed value of hr = 6.2 x 10““^ is good and is an indication of the accuracy of our approximations 
in this region. 

For A = 1, where there is slip on each interface, we hnd that for G < 1, which is always satished, there 
are no real roots for Knioc so that v{x = 0) < 0 for all h. In other words, in this case there is no how 
reversal adjacent to the solid, in agreement with the streamlines shown in Figure [^,c. The accuracy of the 
lubrication approximation for this how is conhrmed in Figure [TO} where the computed velocity tangential to the 
gas’ boundaries Ut is plotted against the vertical velocity predicted by ( [l4| . Note that here, we have used the 
computed value of V at the free-surface, to enable a more accurate assessment of the lubrication approximation, 
rather than the estimate V = 0.85 used elsewhere. 
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Figure 10: Comparison of the computed gas velocity tangential to the boundaries of the gas with the vertical 
velocity from the lubrication analysis for the case ot A = 1 and Pg = 0.014. Curves 1,2 are, respectively, com¬ 
puted values along the free-surface and solid surface whilst corresponding curves from the lubrication analysis 
are la in red and 2a in blue. Notably, there is no flow reversal as u < 0 in all cases. 



Figure 11: Normalised pressure gradient Ap/Apo given as a function of Kuh, where curve 1: A = 0 (in red) 
and curve 2: A = 1 (in black). Circles labelled a, b and c are the values computed for the pressure gradient 
at y = 0.3 for the cases shown in Figure [^,b,c and these compare well to the predictions of the lubrication 
analysis. 


In Vandre et al. (20131, it is shown that air entrainment occurs when the capillary forces at the inflection 
point cannot sustain the pressure gradients required to remove gas from the lubricating gas film. Here, the 
pressure gradient Ap required to maintain the flow at the inflection point, where Kniodh = H) = Kuh, is 
given by 


A - 

oy 


= P 


h=H 


d'^v 

dx^ 


h=H 


6 fl [1 2AKnH + V {2KnH + 1)] 

' m [I + ARnnil -b A) -b UAKnjj]' 


(15) 


As Kuh —>■ 0, so that one approaches no-slip on the solid surface, we have Ap —)• Apo = —6/2(1 -I- V)/H"^ and 
in Figure 11 the pressure gradient in (15), normalised by Apo, is given as a function of Knn for the cases of 
A = 0,1 and shown to agree well with the values from our computations. Notably, curves 1, 2 for A = 0,1 show 
that for large Knu the normalised pressure gradient for A = 0 asymptotes towards a non-zero value whilst for 
A = 1 it tends to zero. Determining Ap/Apo as Kuh —>■ oo from ( |l5| confirms that for A = 0 the limit is finite 
at V/{2{\ + V)) = 0.23 whilst for A = 1 the limit is zero. 

These findings help us to understand the results shown in Figure For A = 1, as the ambient pressure 
is reduced [Pg 0), and hence Maxwell slip is increased {Knn —>■ oo), the pressure gradient required to 
pump gas from the contact line region vanishes, so that the gas offers no resistance to contact line motion and, 
consequently, the maximum speed of wetting appears to become unbounded Cuc —>■ oo. In contrast, for A = 0, 
the lack of slip at the free-surface means that gas is always being driven into the contact line region by the 
motion of the liquid so that however much Pg is reduced {Knu increased) a pressure gradient in the gas is 
always required and Cac approaches a finite value based on this. 
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Figure 12: Deviation of the gas pressure Pg from its far field value of pf as a function of distance from the 
contact line y — Pc with curves 1 (in black) and 2 (in blue) corresponding to, respectively, cases considered in 
Figure and Figure [^. The pressure along the gas-solid interface is given by curves la, 2a whilst 16,26 are 
values from the gas-liquid interface. 


4.5 Evaluation of the Incompressibility Assumption 


The analytic expression derived for the pressure gradient in the thin film (151 can be used to estimate the 


validity of the assumptions we have made about the gas flow. In particular, although the flow is clearly low 


Mach number, compressibility effects in thin films can still be significant, see §4.5 of Gad-el-Hak (20061. This 


occurs when the large pressure gradients required to pump gas out of the thin film lead to reductions in pressure 
that are comparable to the ambient pressure in the gas. As it is the gas flow at the inflection point on the 
free-surface that is most important for air entrainment, we will estimate whether the flow in the gas there is 
indeed incompressible. 

The pressure gradient obtained in (15) takes a maximum value of Apo = — 12/i/iF^, as 14 < 1, so that the 
maximum pressure change along a film of length T/ is 12flLf/H^. Simulations show that the film length is 
never larger than the capillary length so that Lf = 1 is an upper bound. Comparing the pressure change to 
the ambient pressure Pg which (dimensionlessly) is Pg = P*L* /[pPU*) gives us that compressibility can be 
neglected if Apo ^ Pg which requires 


H^Ht = 


3.5^p*U*/iP;L*), 


(16) 


where Ht is the (dimensionless) transition height below which compressible effects can no longer be neglected. 
This is consistent with a similar condition derived for drop impact phenomena in Mani et al. (2010). 

Taking U* = 1 m s“^ gives Ht = 1.2 x 10“^ for air at atmospheric pressure whilst reducing the pressure 
by a factor of one hundred gives Ht = 1.2 x 10“^. When looking at the values of H obtained at various Ca 


in Figure]^ it appears that compressibility may indeed be important. However, the estimate obtained in (16) 
overpredicts Ht due to (a) neglecting reductions in pressure gradients associated with slip at the interfaces 
and (b) approximating the film as being a channel of length L and constant height H whereas the film height 
actually increases by orders of magnitude as one moves away from the contact line, see Figure Therefore, 
before abandoning incompressibility, a more accurate evaluation of this assumption is considered in which the 
pressure changes along the film obtained from our computations are used instead of Apo from 

In Figure 12 the variation in pg from its far-field value p/ is plotted as a function of distance from the 
contact line y — Pc along both the gas-solid and gas-liquid interfaces for the cases shown in Figure [^,c, i.e. for 
Pg = 0.14, 0.014 with A = 1. As one would expect in a lubrication flow, the values along the two interfaces 
are close and are graphically indistinguishable along curve 2. It can be seen that the maximum change in pg in 
the thin film region for the cases considered is 17 and 3, respectively. For the liquid associated with the base 
state, the substrate speed at Ca = 0.4 is U* = 0.16 m s“^ so that for Pg = 0.014, where the gas is most likely 
to be compressible, we have Pg = 260. Therefore, in these cases the pressure change along the film Pg — Pf is is 
substantially smaller than the ambient pressure Pg with {pg — pf)/Pg < 0.012. 

Notably, although the lubrication analysis accurately predicts the pressure gradient in the film, see Figure [TT| 
using the analytic value at the inflection point to estimate the pressure change along the entire length of the 
thin film is inaccurate. This is due to the widening of the film as one moves away from the contact line. A more 
accurate estimate could have been achieved by integrating the analytically calculated pressure gradient over the 
computationally obtained film profile, but a simpler way is to use the computed values of Pg. Doing this shows 
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Figure 13: Dependence of the critical capillary number Ca^ on the viscosity ratio fi, with all other parameters 
fixed at their base values, with curve 1: Pg = I (black), 2: Pg = 0.1 (red) and 3: Pg = 0.01 (blue). 

that the assumption of incompressibility is an accurate one for the base cases considered. Further simulations 
confirm that even in the most extreme cases considered here, the incompressibility assumption remains a good 
one. 


5 Parametric Study of the System 

Having fully analysed the base state, the role of the system’s parameters is now established by perturbing their 
values about the base ones. Henceforth, we only consider the most relevant case of H = 1. 


5.1 Influence of Viscosity Ratio 

The viscosity ratio fl has long been recognised as an important parameter in coating flows as it is a measure 
of how much resistance the receding gas phase can produce. Working with a dimensionless system allows us to 
isolate the effect of ft on Cuc, at different ambient pressures, whilst keeping all other parameters fixed at their 
base state. The range considered is chosen to cover all values of p, obtained for the liquids in |Benkreira fc Khan| 
(20081 used to coat a solid in air, giving 10“^ < p < 10“^. 

In Figure 13 one can see that lower values of the viscosity ratio p result in higher Cqc across all values of 


Pg. What was less expected, is that for Pg = 0.01 (curve 3) it appears there is a value of p below which Cuc 
increases apparently without bound. This occurs at p = 3.3 x 10““^ which is slightly less than the base case’s 
viscosity ratio {po = 3.6 x 10“"*) where the critical gas pressure Pg c = 9 x 10“^. Therefore, at smaller p the 
critical value Pg c increases, i.e. less pressure reduction is required to reach the critical value. 

The effect of p on Pg c can be understood by looking back to the lubrication analysis in the previous section 
and, in particular, the expression for the pressure gradient in the gas (151 which is proportional to p (in 


dimensional terms it would be proportional to fj.*). Therefore, whilst the degree of Maxwell slip controls the 
boundary conditions to the gas flow, the bulk flow in the gas is characterised by the viscosity ratio p, with 
larger values decreasing Cuc as more effort is required to remove gas from the contact line region. 


5.2 Effect of Ohnesorge Number 

The effects of inertia are usually assumed to have a negligible influence on the dynamics of air entrainment, with 
Stokes flow, corresponding to Oh —>■ oo, often considered. In full computations of coating flows in Vandre et al. 


(20131 it was shown that (a) inertial effects in the gas have a negligible influence, reaffirming our assumption 
that p = 0 can be taken without loss of generality, and (b) inertial forces in the liquid do not alter Coc until 
they are relatively large. The results in Figurej^at atmospheric pressure agree with these previous conclusions, 
with a 10^ reduction in Oh only increasing value of Cac by about 10%. Note that Cac at Oh = 0.1 corresponds 
to Re = 52. 

For the lowest ambient pressure Pg = 0.01 an entirely different behaviour is observed and there is a critical 
value of Oh where Cac appears to increase without bound. In contrast to the influences of Pg and p which 
could be rationalised by considering the gas flow, this effect is generated by changes in the liquid’s dynamics. 


In particular, in Vandre et al. (2013) it was shown that when inertial effects become more prominent, velocity 


gradients in the liquid are localised to a thin boundary layer near the moving solid surface and, as a result. 
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Figure 14: The dependence of the critical capillary number Coc on the Ohnesorge number Oh with all other 
parameters fixed at their base values, with curve 1: Pg = I (black), 2: Pg = 0.1 (red) and 3: Pg = 0.01 (blue). 



Figure 15: The dependence of the critical capillary number Coc on the equilibrium contact angle 9e with all 
other parameters fixed at their base values, with curve 1: Pg = 1 (black), 2: Pg = 0.1 (red) and 3: Pg = 0.01 
(blue). 


the free-surface becomes less deformed so that Coc increases. Our results show that when combined with small 
values of Pg this mechanism can prevent air entrainment at any Ca. 


5.3 Role of Substrate Wettability 

Experimentally, the wettability of the substrate is known to have an influence on the point of air entrainment 


both in coating flows and in impact problems, such as the impact of solid spheres on liquid baths (Duez et al. 


20071. In Figure 15 our results obtained at constant contact angles 6 d = 0, 

is. Curves are from 5° < 6 , 


show that the more wettable the 
to avoid computational difficulties associated with 


substrate, the higher Coc is. Curves are from 5° < 9e < 175° 
extremely small angles in either phase, but, despite this, the trends are still clear enough. Notably, at reduced 
pressures variations in Ca^ with 9f, are far more dramatic: changing 9^ from 30° to 90° reduces Ca^ by 0.32 at 
Pg = 0.01 and compared to 0.1 at Pg = 1. 

The model used here predicts that Coc monotonically decreases as 9e increases, which is what one may 
intuitively expect. However, in Blake & De Coninck (2002) the optimal value of 9e which maximises Cuc is 
not at zero and can even occurs on hydrophobic substrates 9^ > 90°. It is likely that such effects can only be 
captured by using a dynamic contact angle model and so this aspect lies beyond the scope of the present paper. 
However, further experiments at both atmospheric and reduced pressure for solids of different wettabilities 
would be useful shed further light on these issues. 
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5.4 Different Gases 


Previous studies in both dip-coating (Benkreira & Ikin 2010) and drop impact (Xu et al.. 2005) have considered 


the effect of using different gases on wetting failure. The gases used tend to have a similar viscosity and density 
at atmospheric pressure to air but possess different molecular weights and hence mean free paths. The fluid 
flow considered here is incompressible, so that changes in the speed of sound in the gas, and hence the Mach 
number, do not concern us. However, changes in the mean free path at a given ambient pressure will have an 
effect on our results. This effect can be accounted for in the results presented, without recomputing everything, 
by simply rescaling the ambient gas pressure Pg. For example, if a new gas considered has a mean free path 
£atm,n which is c times larger than that of air at atmospheric pressure, so that iatm,n = then from (|^ 

the same Kn is obtained if we use a rescaled pressure Pg^n in the new gas which satisfies Pg^n = cPg. 

In this way, all the results of the previous section can be used for any gas, rather than just air. For example. 


if helium is used instead of air, approximately 


atm.He 


— 0-c„ 


atm,air’ SO that Pg,He = ^Pg- Then, for the base case 
considered in 14.1 the critical pressure was found to be Pg = 9 x 10“^ for air so that for helium this would be 

Therefore, as seen experimentally, built into the model is the fact that gases with 


increased to Pg = 2.7x10 ^ 


larger mean free paths require less pressure reduction in order to induce signihcant non-equilibrium gas effects 
that increase Car- 


5.5 Summary of the Parametric Study 

The effects of non-equilibrium gas dynamics on air entrainment phenomena have been clarified and the role 
of the various dimensionless parameters on Cgc has been characterised by perturbing them about a base case. 
This has allowed us to isolate the effects of parameters such as p, which cannot easily be independently varied 
experimentally. The disadvantage is, of course, that the results cannot easily be used to find, say, the effect of 
fj.* on coating speeds as this variable comes into Ca, Oh and p which have thus far been varied independently. 


Therefore, in the next section we compare our results directly to experiments in Benkreira & Khan (2008) and, 
in doing so, now step into a dimensional setting. 


6 Comparison to Experimental Data: Effect of Viscosity and Gas 
Pressure 


Consider now whether the new model is able to account for the influence of gas pressure P* on maximum 
coating speeds U* observed in Benkreira & Khan (2008) for different viscosity liquids. All parameters remain 


unchanged from ^except those which depend on viscosity, namely Oh, Ca and p which will be based on values 
used in the experiments of /x* = 20, 50, 100, 200 mPa s. 

In Figure 16 a direct comparison between the computations and the experimental results in Benkreira & 
is shown. 


Khan (2008) is shown. Notably, both predict that a significant reduction in P* from its atmospheric value 
(Pg,atm = 1^0 kPa) is required to induce increases in U*, but that once this has been achieved, the subsequent 
enhancements in U* are substantial. At the lowest viscosity, the agreement is relatively good throughout; 
however, as is most clear from Figure [T6|3, one can see that at the higher viscosities, the effect of gas pressure 
yields enhanced coating speeds earlier in the experiments than in the computations. For example, for the highest 
viscosity liquid, at P* = 12 kPa the experimental data (crosses) shows that U* has increased to 0.18 m s“^ 
from its atmospheric value of 0.11 m s“^ whilst the computations, shown by curve 4, only give this degree of 
enhancement once P* = 2 kPa. 

Promisingly, as in the experimental results, curves 1 and 2 cross at low atmospheric pressure meaning that 
a higher viscosity liquid can be coated at faster speeds than a lower one, in contrast to what is observed at 
atmospheric pressure. Other crossovers may have been observed if higher capillary numbers could have been 
obtained, as suggested by the path of curve 3. However, the position of the crossovers do not match the 
experimental data; whilst the computed results all appear to crossover around P* ~ 1 kPa, the experimental 
results do so at higher values of P* with the two highest viscosity solutions crossing over at P* ~ 10 kPa. 

In summary, whilst all of the qualitative features of the experimental results are recovered in our computa¬ 
tions, quantitative agreement is not good over the entire range of data. In particular, less pressure reduction 
is required in the experiments than in the computations for there to be a noticeable effect on the maximum 
wetting speed. There are a number of possible reasons for this discrepancy, some of which we will now consider. 


Roughness of the Solid: In Benkreira & Khan (2008) the “average peak-to-valley height roughness” for 
the experimental results presented is 300 nm. Although this is ‘smooth’ compared to the other solids used 
in this work, it is still quite rough compared to the solids used in some other coating flow studies. This 
effect cannot be ignored as experiments in [Clarke (2002) have demonstrated that the roughness of a solid 
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Benkreira & Khan 


Figure 16: A comparison of our computations for the critical wetting speed U* as a function of gas pressure 

pOMj). Curves are computed for liquids of different 
ck), 

(in green) and the corresponding experimental data are given by, respectively, triangles, diamonds, squares and 


Pg to experimental results obtained in 
viscosities with 1: 20 mPa s (in red), 2: 50 mPa s (in black), 3: 100 mPa s (in blue) and 4: 200 mPa s 


crosses. 
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can influence the maximum speed of wetting with, in many cases, faster coating speeds obtained with 
rougher surfaces. 

From a modelling perspective, the level of roughness In the experiments is large when compared to the slip 
length at the liquid-solid interface I* = 10 nm. However, simply increasing the value of I* In an attempt 
to mimic the effect of the roughness does not capture the experimental data either. It is likely that the 
strongest effect of the roughness will be on the actual contact angle, but a theoretical description of this 
mechanism at high contact-line speeds remains an open problem and is beyond the scope of this work 


Dynamics of the Contact Angle: Despite assuming that the contact angle Is fixed at its equilibrium value, 
without fitting the slip length in any way we still obtain reasonable agreement for the values of U* at 
atmospheric pressure. This could be a coincidence and dynamic contact angle effects may Indeed be 
important for quantitatively capturing this class of coating flows; however, such effects will not Improve 
the agreement of our current model with experiments. This is because allowing the contact angle to vary 
will Increase the angle from its equilibrium value so that, as shown by the findings of |5.3[ U* will decrease 
for a given pressure. Therefore, although the dynamics of the contact angle will alter our results, it will 
not explain why Coc increases so quickly with reductions in Pg and is thus not our main concern here. 


Model for Non-Equilibrium Effects: As shown in §4.3[ the Knudsen number characterising the gas film 
flow indicated that its dynamics is usually outside the ‘slip regime’ as air entrainment is approached and 
more complex non-equilibrium gas dynamics should be incorporated. As noted in §4.3[ incorporation of 
these effects is likely to result in the same qualitative trends but enhanced quantitative ones. In particular, 
it is likely that less gas will be driven into the contact line region at a given pressure so that there is an 
enhancement in Cqc at higher ambient pressures, as observed experimentally. Furthermore, this effect 
will be particularly important at smaller /2, where most of the current discrepancies between theory and 


experiment appear, as the smaller values of H obtained there (Vandre et al. 
importance of non-equilibrium effects 

Models capable of predicting these effects will be discussed in 


2013) will increase the 


The most likely cause of the observed discrepancy seems to be our inadequate description of the gas dynamics, 
and this aspect will be discussed further in ^ Before doing so, we consider what consequences our results have 
for the drop impact phenomenon. 


7 Non-Equilibrium Gas Effects in Drop Impact Phenomena 


It has previously been suggested ( [Rein &: Delplanque |2008[ ), that certain classes of splash observed in drop 
impact phenomena are triggered by the same mechanisms that cause wetting failure in coating flows. This 
interpretation is supported by experiments such as those in Driscoll & Nagel (2011) who find “no trapped air 


beneath the spreading drop outside the small central bubble; there is no significant air film beneath the drop at 
the time of thin-sheet ejection. This suggests that, rather than an underlying air layer, gas flow at the edge of the 
spreading drop is responsible for destabilizing the liquid.” This mechanism is shown by Figure [^. Therefore, 
we now consider how our results can be related to the more complex unsteady drop impact phenomenon. 


The publication which opened up this field of research was Xu et al. (20051. This work focuses on a threshold 


ambient pressure which must be achieved in order to suppress splashing for a particular liquld-solld-gas 
configuration. What is particularly interesting Is that for a given impact speed Vq, the threshold pressure 
strongly depends on the gas, with air, helium, krypton and SFg all used. The data Is shown in Figurewhich 
corresponds in Xu et al. (2005) to the inset of their Figure 2b. In Xu et al. (2005), compressibility effects in the 


gas are assumed to be important when deriving a splash threshold, so that changes In the molecular mass of 
the gases to*, which alter the speed of sound in the gas, become important and result in a collapse of the data 

when P^ is scaled by where to*,^ Is the molecular mass of air. 

An alternative explanation is that compressibility effects are negligible and that Increased Maxwell slip 
is responsible for splash suppression at reduced pressure. Given that the mean free path in the gas at the 
threshold point is inversely proportional to both the gas pressure and the square root of the molecular mass 
to achieve a given increases in to* require decreases in P^. In other words, gases with 


oc (Py 

higher molecular mass will need lower pressures to suppress splashing, as seen in IXu et al.\ (|2005|). 


The dimensionless number characterising the level of Maxwell slip at a given threshold pressure in the gas 
is the Knudsen number Kut = P^/R*, where R* = 1.7 mm is the drop radius. Using the properties of air, the 
expression for Krir required is 


Kut = 


Kn. 


air,atm 


where 


Pt = 


p* 

J /t+'i 


and 


TO„ 


(17) 
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Figure 17: Critical gas pressure required to suppress the splash of a drop with impact speed Vq with different 
gases helium (open circles in green), air (filled circles in blue), krypton (red crosses) and SFg (triangles in black). 
Obtained from Xu et al. (2005). 
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Figure 18: Critical Knudsen number Ktit, scaled by its value in air at atmospheric pressure, required to suppress 
the splash of a drop whose capillary number based impact speed is Cay*- The collapse of data from Figure 17 
onto a single curve is apparent. 


where the coefficient of proportionality Knair,atm = 4.1 x 10“® is the Knudsen number of air at atmospheric 
pressure and the molecular weights in Daltons of = 4, = 29, = 83.8 and rrigp^ = 146 can be 

taken from Xu et al. (20051. As for the coating flow analysis, the value of the Knudsen number has to be based 
on the global characteristics of the flow and would be larger if it could have been derived from the appropriate 
length scales near the moving contact line. 

Ideally, the collapse of data would use Knx and a capillary number based on the contact line speed U* at 
the point of air entrainment, as considered for the coating flow analysis. However, such measurements were not 
made in the experimental work so we shall instead base the capillary number Cay* on the impact speed Kg*, 
an assumption which will be discussed in 

the raw data from Figure 


In Figure 


18 


17 


is shown to collapse well onto axes of Cay* vs Kut- This is 


encouraging, and suggests that non-equilibrium effects in the gas may be important for drop impact phenomena. 
It is particularly reassuring to see that the capillary numbers obtained are within the range one may expect for 
air entrainment phenomena to be occurring, although the capillary number based on the contact line speed at 
the entrainment point may be much larger. Ideally, changes in the viscosity of the liquid would also be included 
into the data collapse, as achieved in Xu et al. (20051, but this is not possible here as such variations cannot 
simply be absorbed into the capillary number, but will also alter many other dominant parameters such as the 
viscosity ratio. 

Although the data collapse looks impressive, a number of outstanding questions remain. First, as also seen 
in m experimentally only relatively small changes in atmospheric pressure are required to suppress splashing, 
with a factor of ten decrease sufficient to prevent it in all cases considered here (Figure [T^, and it remains to be 
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seen whether the Maxwell-slip model is able to predict such a large effect at these reduced pressures. Second, 
in this framework it is not obvious why the drop spreads for some distance before splashing, when the largest 
contact-line speeds are obtained earliest in the impact process. These two issues will be addressed in ^ and 
will require further experimental and theoretical work before they can be fully resolved. 


8 Discussion 

It has been shown that non-equilibrium effects in a gas, influencing the flow through Maxwell slip at boundaries, 
can substantially increase the critical speed at which air entrainment occurs. Notably, neglecting to use slip 
at the gas-liquid boundary produces qualitatively different results from those obtained when slip is allowed, 
predicting a finite, rather than unbounded, maximum speed of wetting as the pressure is reduced. Importantly, 
the problem has been formulated locally, so that the model is valid for any shape of the gas domain and remains 
accurate both far away from the gas film, where non-equilibrium gas effects will have a negligible influence on 
the overall flow, as well as very close to the contact line where regions of high curvature invalidate the lubrication 
approximation. 

proposed was shown to reproduce many of the features of the experimental results in [Benkreii^ 
) but failed to accurately predict the ambient pressures at which coating speeds begin to rise 
for the higher viscosity liquids. These observations open up a number of potential areas for experimental and 
theoretical research which will now be considered. 

8.1 Extension of the Non-Equilibrium Gas Model 

General research into non-equilibrium gas dynamics has identified different ‘regimes’ of flow which indicate the 
modelling approach that should be used. Applying these to the flow in the gas film, characterised by Kriu, 
gives that : 


The model 


& Khan (2008 


Conventional regime {Kuh < 10 ^): The Navier-Stokes equations with no-slip at solid surfaces provide 
an accurate representation of the gas flow. 

Slip regime (10“^ < Kuh < 10“^): Non-equilibrium effects become significant close to the boundary so 
that the no slip condition must be relaxed to allow for Maxwell slip. 

Transition regime (10“^ < Knn < 10): Knudsen layers occupy significant proportions of the domain and 
kinetic theory is required to provide reliable predictions. 

Free molecular flow regime {Kuh > 10): Molecular collisions become rare so that the collision operator 
in the Boltzmann equation can be neglected. 


Clearly, identifying which regime is important for the dynamic wetting problems of interest will be key to 
deciding which methods should be implemented. 

The Knudsen numbers recovered for the gas film flow show that in many cases the model proposed in 
0 has been used outside its limits of applicability. In particular, once in the transition regime the effect of 
Knudsen layers can no longer be accurately captured by a first-order slip condition, as assumed here. An entire 
hierarchy of approaches of varying complexity have been proposed in the published literature to overcome this 
inadequacy, and which of these is optimum for the class of problem considered here is, at present, unclear. A 
sensible first step would be to remain within a continuum framework and to incorporate more advanced non 
equilibrium gas effects either, in the simplest case, by using a second-order slip condition, such as that proposed 


ir 

Hadjiconstantinou 

(2003 

), or, at the next level of difficulty, by modifying the bulk constitutive relations also 

(I 

jockerby et al. 

2005 a 

). Research in this direction is currently underway. 


The regions of applicability of the aforementioned approaches are often unclear, and the only way to test these 
are to compare with solutions obtained from the Boltzmann equation. This can now be routinely achieved using 
a variety of methods, with the direct simulation monte carlo (DSMC) approach popular (Bird 19941. However, 
in this setting it is not clear how the liquid-gas free-surface should be treated as the Boltzmann equation 
cannot be applied across the boundary and into the liquid phase. Whether some plausible assumptions about 
the liquid’s free-surface properties could be made or whether DSMC in the gas could somehow be coupled to 
molecular dynamics simulations (Koplik & Banavar, 19951 in the bulk of the liquid remains to be seen. Progress 
in the latter area would shed light on not only the bulk flow of the gas, but also its interaction with boundaries 
formed by liquids. 
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8.2 Alternative Flow Configurations 


It is of interest to look at the effects of gas pressure on coating speeds in a range of different setups, both 
experimentally and theoretically. In particular, the effects of confinement, where the contact line is isolated by 
a nearby boundary, have not been considered and neither have axi-symmetric geometries such as those used in 
fibre coating. 


The fibre coating configuration is particularly interesting as experiments in Simpkins & Kuck (20031 show 
that extremely high capillary numbers can be achieved before air entrainment occurs and, at present, the 
precise cause for this remains to be seen. Some works claim that the same mechanism has been observed in 


hydrodynamic assist (Blake et al.. 2004) whilst others have proposed that air is actually entrained into the liquid 


but then dissolved due to the high pressures there (Jacqmin 20021. Incorporating the models proposed into the 


computational framework developed here should allow one to distinguish between these rival interpretations. 


8.3 Development of Analytic Theory 


The theory of |Cox| ( |1986[ ) has widely been used to calculate an apparent contact angle from its actual value 
in both single and two-phase dynamic wetting phenomena. Although this theory is based on a small capillary 
number asymptotics, so that it lacks accuracy for many coating flows (Vandre et at, 2013), it nevertheless has 


been a useful tool for interpreting experiments and predicts a possible mechanism by which air entrainment can 
occur; namely, that once the apparent angle reaches 180°, no more solutions can be obtained. However, in its 
present state this theory cannot capture the phenomena investigated in this thesis. 

Cox’s theory contains an ‘inner region’ where additional effects come into play, such as slip, that allows the 
moving contact line problem to be circumvented. When one has the same slip length in each phase, then the 
theory is valid, but for the problems we have considered the slip in the gas phase is many times larger than that 
in the liquid. In particular, the slip length in the gas phase is likely to become comparable with the scale of the 
‘intermediate region’, where the bending of the free-surface occurs, so that the entire asymptotic structure of 
the problem becomes unclear. It would be of interest to see if Cox’s approach can be generalised for this more 
complex class of flows. 


8.4 Parametric Experimental Study of Dip Coating at Reduced Pressures 


The results of Benkreira fc Khalil ( |2008[ ) have allowed us to compare our computational results to experimental 
data. However, these are, at present, the only results in the published literature considering the effect of gas 
pressure on coating processes. Therefore, additional experimental studies would be useful in order to develop 
the theory of this class of problems. Such an investigation may, for example, use a wider variety of liquids 
and consider solids with a range of different wettabilities and/or roughnesses. Furthermore, considering how 
accelerating or decelerating the plate alters the flow may give some initial insight into more complex unsteady 
flows such as drop impact. 

If the aforementioned analyses could be combined with modern experimental techniques, such as interfer¬ 
ometry measurements, then some details of the gas film’s profile could potentially be recovered, as has been 
achieved in drop impact phenomena. These results would provide an excellent benchmark for any potential 
theories to be compared to. 


8.5 Drop Impact and Spreading 

Drop impact phenomena are notoriously difficult to capture computationally due to the inherently multiscale 
nature of a problem which involves unsteady large deformation of the flow domain. Although this has been 
achieved in the case of a single-phase flow, e.g. in Sprittles & Shikhmurzaev (20I2&), at present there are no 


codes which have captured all physical scales in (a) the dynamics of impact, where a thin gas film cushions the 
falling drop; (b) the topological change, when the drop first touches the solid; and (c) the spreading process, 
where the slip length needs to be resolved. Building a code capable of achieving this goal is highly sought after. 

A question arising from our analysis in 0 is why the drop spreads for some distance before splashing, when 
its highest speeds are obtained early on in the impact process. The answer to this question is likely to involve 
the drop’s shape and its flow field, so that one cannot simply consider an isolated moving contact line, as we 
have done, but must also consider the outer influence of the flow field on the point at which air entrainment 
occurs. 

The situation with drop impact could be similar to that of curtain coating where the point at which air 


entrainment occurs is dependent not only on the contact line speed, but also on the flow geometry (Blake et al. 
1994 Clarke & Stattersfield 2006). In particular, when the contact line is below the falling liquid curtain 


it appears to be ‘assisted’ by the downward pressure of the curtain (]Blake et al. 2004), which prevents air 


entrainment, whilst if a heel forms in the free-surface the contact line becomes isolated from the bulk flow and 
there is no assist. Similarly, in drop impact, during the initial stages of spreading the contact line is located 
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below the drop, and thus wetting is assisted by the downward pressure of the bulk flow, until a lamella is ejected, 
analogous to a heel in curtain coating. At this point, the contact line is isolated and air entrainment between 
the lamella and the substrate is able to lift off a sheet of liquid which creates a splash. 

In order to determine if the aforementioned mechanism is the correct one for drop impact phenomena, it 
may still be useful to use the simpler coating flow setups. One possibility is that the influence of an impulsively 
started or accelerating plate could be considered in a manner to mimic the impact process. Another avenue of 
enquiry would be to alter flow domain to mimic the flux of mass arriving at the contact line in drop impact 
phenomena. 

From an experimental perspective, it would be extremely useful to know the contact line speed at the point 
of air entrainment, so that the capillary numbers used in ij^can be based on this, rather than the impact speed. 
In particular, if air entrainment occurs after the drop’s contact line has moved a radius r*, and the contact line 


moves with a square root in time dependency, as shown experimentally (Riboux & Gordillo, 20141 for the initial 
stages of motion, so that r*/R* = A^V^t* jR* and U* jV^ = A/2^R* /{VQt*), where A is a constant, we find 
that U*/Vq = {A^/2)R*/r*. Therefore, Cay* is only a fixed multiple of Ca if air entrainment always occurs 
at the same radius r*, which seems unlikely to be the case across all parameter space. Therefore, we cannot be 
sure that Cuc oc Cay* and more experimental data is indeed required. 


Appendix: Benchmark Calculations 


From a theoretical perspective, both extended lubrication theory (Chan et al. 20131 and full computations 


(Vandre et al.. 2012 20131 have identified a critical capillary number past which no steady two-dimensional 


solutions exist, and it has been assumed that this corresponds to the Cuc observed in experiments. As initially 


observed for the plate withdrawal problem, i.e. the dewetting scenario (Snoeijer et al. 20071, the turning point 


at Cqc in the map of elevation Ay (or maximum apparent contact angle 9m) versus Ca, separates stable and 


unstable branches of steady-state solutions (Figure 191. Consequently, when considering the entire solution 


space, for certain Ca multiple values of elevation Ay exist, and no function Ay = Ay{Ca) exists. Notably 
though, along all sections of the solution path we have Ca < Cac, so that Cac is indeed the maximum speed of 
wetting. 

In order to determine Cac for a given system, a number of possible solution techniques exist. The most 
obvious approach is to find a solution at a relatively small Ca, then increase its value, using the previous 
computed solution as an initial guess for the next Ca, and repeat until no more steady state solutions exist. 
This is the simplest method and reliably finds Cac, but it does not allow us to extract the solution path past 
the critical point where Ca < Cac, as noted in Vandre et al. (2012). A different approach, which allows the 


entire solution path to be captured, is to determine Ca, i.e. treat it as an unknown in the problem, for a fixed 
elevation of the contact line Ay. This approach works better in practice, as for every Ay, there is a unique Ca, 


i.e. there is a function Ca = Ca{Ay). Consequently, curves like those in Figure 19 can be generated. 


In principle, as much of the unstable branch as required can be captured, but for Ay Aycrit, where Aycrit 
is the elevation when Ca = Cac, the deformation of the free-surface becomes large and this process becomes 
computationally intensive, particularly at small /2. Consequently, simulations were run up to the point in the 
solution path where the mesh failed, usually due to the finite elements becoming too distorted, and stopped 
there. As one can see from Figure [T^ all simulations were extended well past Ay = Aycrit and it was only the 
smallest viscosity ratio (curve 4) where the code finished prematurely, i.e. before Ay = —3. 


Calculations and Comparison with Vandre et al. (2013) 


In order to compare our computations to previous work in Vandre et al. (2013), consider Stokes flow Oh~^ = 0; 
with the slip-lengths of the liquid-solid and gas-solid interfaces equal at Is = Kn = 10“"^; no-slip across the 
liquid-gas interface, so that A = 0; and a constant contact angle of 9d = 90°. Using these parameters, the critical 
capillary number Cac and maximum apparent contact angle at this point 9m, which occurs at the inflection 
point on the free-surface (see F igure [^, are co mputed for different viscosity ratios fi = 10“^, 10“^, 10“^, 10“^ 
and compared to the results in Vandre et al. (2013). 


In Figure [T9l the dependence of Ay and 9m on Ca are computed. From the curves one can clearly see the 
existence of a critical capillary number Cac past which (Ca > Cac) no solutions exist. On further decreasing 
Ay, so that the contact line moves further down the plate, solutions are recovered for which Ca < Cac, and 
the capillary number begins to oscillate around a capillary number Ca* < Cac as also observed in theoretical 
works on dewetting phenemena, where these solutions have been shown to correspond to a series of saddle-node 


2006). 


bifurcations, see Figures 4 and 5 of Snoeijer et al. (2007), which are related to transient states (Snoeijer et al. 


Similar oscillations are observed for 9m where, notably, along the unstable branch the maximum apparent 
angle can increase beyond 180°. This can be seen in Figure!^ where between curves 3 and 4 a hump develops 
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Figure 19: Dependence of the elevation Aj/ and maximum apparent contact angle 9m on the capillary number 
Ca for p, = 1 : 10“^, 2 : 10“^, 3 : 10“^ ,4 : 10“^ and other parameters listed in the text. Dashed lines indicate 
the critical capillary number Ca = Cuc past which no steady state solutions exist. 
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Figure 20: Free-surface profiles for the case of /x = 10 ^ at different elevations 1: Ay = —0.5, 2: Ay = —1, 
3: Ay = —1.5, 4: Ay = —2, 5: Ay = —2.5 and 6: Ay = —3. 




CcLq 


9M{Cac) 






Computed 

Vandre et al. (2013) 

Computed 

Vandre et al. (2013) 

10-^ 

0.09 

0.09 

161° 

159° 

10-2 

0.21 

0.22 

165° 

166° 

10-3 

0.39 

0.41 

171° 

172° 

10-4 

0.65 

0.65 

176° 

176° 


Table 2: Comparison of the critical capillary numbers Coc and free-surface angle at the inflection 
obtained from our code with those obtained in Vandre et al. (20131 across a range of viscosity ratios 


point 9m 

fl. 


in the free-surface profile. In this region, the gradient of the free-surface is negative so that the angle which this 
portion makes with the solid, i.e. the apparent angle, gives a maximum of 9m > 180°. 

A comparison between the critical capillary numbers computed in Vandre et al. (2013) and those observed 
here is not likely to be exact, as different flow geometries have been used. However, in cases where Cqc is not 
too large and the effect of the flow geometry is relatively weak we may expect our values to be comparable to 
theirs. In Tablej^this is shown to be the case. Computed values of both Cac and 9m at this capillary number, 
which is a very sensitive measure, agree well across a large range of viscosity ratios. It is not our intention to 
analyse this data in detail, as this has already been achieved in Vandre et al. (2013), but this table will be useful 
as a benchmark for future computational work in this area. 
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